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The  analyses  of  reservoir  flushing  and  salinity  intrusion  events  for  the 
proposed  Bushy  Park  Entrance  Canal  relocation  were  performed  for  the  U3  Army 
Engineer  District,  Charleston  (SAC).  Studies  were  performed  in  the  Hydraulics 
Laboratory  (HL)  of  the  US  Army  Engineer  Waterways  Experiment  Station  (WES) 
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analytical  salinity  model  used  in  this  study.  Dr.  Kurt  Getsinger,  Environ¬ 
mental  Laboratory,  WES,  conducted  tests  on  dye  uptake  by  water  hyacinths. 
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report. 
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CONVERSION  FACTORS,  NON-SI  TO  SI  (METRIC) 
UNITS  OF  MEASUREMENT 


Non-SI  units  of  measurement  used  in  this  report  can  be  converted  to  SI 
(metric)  units  as  follows: 


Multiply _  Bv  To  Obtain 


cubic  feet 

0.02831685 

cubic  metres 

cubic  yards 

0.7645549 

cubic  metres 

feet 

0.3048 

metres 

miles  (US  statute) 

1.6093 

kilometres 
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EFFECTS  OF  A  PROPOSED  BUSHY  -PARK  ENTRANCE  .CANAL 
RELOCATION.  COOPER  RIVER.  SOUTH  CAROLINA 

PART  I :  INTRODUCTION 

BagHerpjaKi 

1.  In  1985,  the  US  Army  Engineer  District,  Charleston,  reversed  a  1942 
flow  diversion  and  rediverted  most  of  the  442-cu  m/sec  Cooper  River  inflow 
back  to  the  Santee  River,  leaving  about  127  cu  m/sec  weekly  average  inflow 
from  the  Pinopolis  Dam.  The  purpose  of  the  rediversion  was  to  reduce  shoaling 
in  Charleston  Harbor,  and  initial  indications  are  that  expected  shoaling 
reductions  will  be  met  or  exceeded.  Since  rediversion,  another  concern  of  the 
Charleston  District  has  become  the  protection  of  the  Bushy  Park  Reservoir  from 
occasional  salinity  intrusion  into  the  entrance  canal.  Bushy  Park  is  a  tidal, 
freshwater  reservoir  located  in  the  upper  reaches  of  the  Cooper  River  (Fig¬ 
ure  1).  About  a  dozen  intrusion  events  have  occurred,  generally  correlated  to 
times  of  tidal  high  water  and  prolonged  zero  discharge  from  Pinopolis,  high 
tidal  amplitudes,  and/or  fluctuations  in  mean  tide  level.  An  alert  system  is 
in  place  to  sense  salinity  intrusions,  and  to  order  supplemental  freshwater 
inflow  releases  by  the  State  of  South  Carolina,  the  operators  of  the  hydro¬ 
plant  at  the  Pinopolis  Dam.  As  a  permanent  structural  solution,  the  Charles¬ 
ton  District  proposed  a  new  facility  to  be  added  to  the  authorized  Cooper 
River  Rediversion  Project  consisting  of  a  canal  to  relocate  the  entrance  up¬ 
stream  to  a  point  above  the  limits  of  salinity  intrusion  (US  Army  Engineer 
District,  Charleston,  1983).  A  concern  in  the  design  of  a  relocated  entrance 
canal  is  its  possible  effects  on  the  flushing  characteristics  of  the 
reservoir. 

2.  The  US  Army  Engineer  Waterways  Experiment  Station  (WES)  has  assisted 
the  Charleston  District  on  a  number  of  studies  concerning  the  Charleston  Har¬ 
bor  and  Cooper  River  in  support  of  the  Cooper  River  rediversion.  WES  first 
became  involved  during  the  mid-1940's  as  consultant  on  the  entrance  to  the 
Bushy  Park  reservoir  at  about  river  mile  43  on  the  Cooper  River.  During  the 
1950's  through  the  1970' s,  two  physical  hydraulic  models  were  used  at  WES  to 
conduct  a  series  of  studies  on  both  Charleston  Harbor  and  Bushy  Park  Reservoir 
conditions  related  to  the  Cooper  River  rediversion  (Benson  1976;  Benson  1977; 
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and  Benson  and  Boland  1977).  Charleston  Harbor  conditions  were  monitored  by 
WES  both  before  and  immediately  after  the  Cooper  River  rediversion  in  1985, 
and  again  during  1987  for  comparison  analysis  (Teeter  1989).  The  studies 
addressed  mixing  and  sediment  flushing  in  the  harbor,  as  well  as  the  inflow, 
tidal,  and  meteorological  effects  on  salinity  conditions  within  the  Bushy  Park 
Reservoir.  A  schematic  numerical  model  was  developed  to  examine  salinity 
intrusion  and  applied  to  predict  harbor-deepening  shoaling  effects  (Teeter  and 
Pankow  1989).  The  US  Geological  Survey  (USGS),  Columbia,  SC,  office  has 
installed  and  maintains  a  system  of  satellite-telemetered  water  level  and 
water  quality  monitors,  providing  real-time  data  for  a  WES-designed  salinity 
alert  system.  Background  salinity  conditions  were  previously  compiled  for  the 
reservoir,  which  the  Corps  is  committed  to  protect  from  salinity  intrusion. 


3.  The  purpose  of  this  study  was  to  evaluate  the  proposed  Bushy  Park 
entrance  canal  relocation  and  corresponding  effects  on  tidal  volume  exchanges, 
flushing  rates,  and  salinity  intrusion.  Dye  tests  were  performed  to  gage 
flushing  rates,  tide  heights,  and  volumes  for  the  present  canal  configuration; 
to  establish  baseline  conditions  for  the  reservoir;  and  to  provide  information 
for  evaluation  of  reservoir  flushing  processes.  The  purpose  of  the  tidal 
model  tests  was  to  evaluate  the  Bushy  Park  Reservoir  entrance  canal  relocation 
and  corresponding  tidal  volume  exchanges  and  to  develop  an  entrance  canal 
geometry  that  would  maintain  tidal  flushing  at  its  current  level.  The  purpose 
of  salinity  analysis  was  to  evaluate  the  effectiveness  of  the  proposed  en¬ 
trance  canal  location  in  protecting  Bushy  Park  Reservoir  from  ocean  salinity 
intrusions . 


4.  The  following  overall  approach  was  used  in  this  study: 

Field  measurements  were  performed  with  USGS  to  obtain  new  data 
(dye  dispersion,  tidal  volumes,  salinity,  etc.)  in  the  existing 
entrance  canal  for  two  tidal  conditions.  This  report  briefly 
describes  prototype  data  collection  activities  and  presents 
selected  field  data  and  an  analysis  of  reservoir  flushing  using 
a  simple  model.  The  model  considered  tidal  exchanges  and 
industrial  withdrawals. 
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b.  Numerical  modeling  was  performed  using  an  established  numerical 
model  to  determine  the  effects  of  existing  and  proposed  entrance 
canal  geometry  (length  and  cross  section)  on  reservoir  flushing 
and  tidal  volumes .  This  report  summarizes  numerical  model  veri¬ 
fication  and  tests  with  alternate  entrance  canal  dimensions  and 
supplemental  withdrawals.  Observed  and  computed  reservoir  tidal 
volume  exchanges  and  water  levels  are  presented. 

c.  Specific  conductance  data  were  compiled  from  historical  records 
to  test  rediversion  salinity  effects  at  the  DuPont  intake  sta¬ 
tion  located  on  the  reservoir  entrance  canal  and  to  extrapolate 
data  from  other  stations  upstream  to  determine  extreme  transient 
conditions  at  the  proposed  entrance  canal  location.  Historical 
data  were  compiled  for  times  of  high  chloride*  readings  in 
Bushy  Park  between  1  January  1986  and  1  August  1989.  Select 
intrusion  events  were  analyzed  using  an  analytic  model  to  pre¬ 
dict  the  maximum  chloride  intrusion  at  the  proposed  entrance 
canal  location. 


The  organization  of  the  next  three  parts  of  this  report  follows  this  outline. 


*  Salinity  is  a  measure  of  the  total  dissolved  salts  in  water,  with  ocean 
water  having  a  salinity  of  about  34  ppt.  Chloride  content  is  a  more  consis¬ 
tent  measure  in  low-salinity  waters  and  is  the  measure  of  preference  for 
Bushy  Park  water  quality  (paragraph  58).  Depending  on  the  context,  both 
measures  are  used  here. 
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PART  II:  PROTOTYPE  SYSTEM  FLUSHING 


Site  Description 

5.  Bushy  Park  Reservoir  is  located  approximately  27  km  north  of 
Charleston,  SC.  Figure  1  shows  the  Cooper  River  and  the  Eushy  Park  Reservoir. 
The  reservoir  was  completed  in  1956  by  closing  the  Back  River  with  an  earthen 
dam  and  dredging  an  intake  canal  to  bring  fresh  water  from  the  Cooper  River. 
Bushy  Park  Reservoir  is  used  primarily  for  industrial  and  municipal  water 
supply. 

6.  The  reservoir  is  about  14.5  km  long  from  the  existing  entrance  to 

the  dam.  The  reservoir  is  surrounded  by  marsh  except  at  the  industrial  site 

between  the  reservoir  and  the  Cooper  River.  The  entrance  canal  intersects  the 
Cooper  River  at  river  mile  (RM)  43  as  measured  from  the  harbor  jetties.  The 
proposed  entrance  canal  would  add  about  3.2  km  to  the  overall  length  of  the 
reservoir,  and  would  situate  the  reservoir  entrance  about  6.4  km  farther  up¬ 
stream  on  the  Cooper  River.  Figure  1  shows  a  schematic  of  the  proposed  canal. 

7.  The  existing  entrance  canal  has  dimensions  of  about  55  m  wide  by 

4.2  m  deep  to  mean  tide  level  (mtl) ,  and  has  been  stable  since  its  construc¬ 
tion.  No  dredging  has  ever  been  required,  nor  are  there  any  indications  of 
shoaling  in  the  canal.  Canal  side  slopes  range  from  IV  on  2H  to  IV  on  4H. 

8.  Tides  are  semidiurnal  with  a  mean  tide  range  of  1.62  m  at  the 
Customs  House  in  Charleston  Harbor.  The  average  tide  range  near  the  entrance 
to  Bushy  Park  Reservoir  on  the  Cooper  River  is  about  0.89  m,  and  just  inside 
the  reservoir  about  0.55  m. 

9.  The  volume  and  surface  area  of  the  reservoir  were  estimated  at 

10  million  cubic  metres  and  2.5  million  square  metres  at  mtl.  For  comparison, 
the  volume  and  surface  area  of  the  existing  entrance  canal  were  estimated  at 
0.8  million  cubic  metres  and  0.2  million  square  metres  at  mtl. 

10.  The  withdrawals  from  Bushy  Park  Reservoir  are  substantial.  Total 
withdrawal  capacity  is  about  31  cu  m/sec.  The  South  Carolina  Electric  and  Gas 
(SCE&G)  withdrawal  alone  is  about  24  cu  m/sec,  about  77  percent  of  the  total. 
Characteristic  times  for  50  percent  replacement  for  the  reservoir  were  roughly 
estimated  based  on  plug-flow  withdrawal  Tw  and  exponential  tidal  replacement 
times  Tt  as  follows: 
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3.6  tidal  cycles 


(1) 


T.-0.5  “ 


Tt  *  -In  (0.5)  ^5  =5  tidal  cycles 


(2) 


where 

V0  -  volume  of  the  reservoir  to  mtl  (without  the  entrance  channel 
-  10  million  cubic  metres) 

Q„  -  withdrawal  rate 

Vp  -  tidal  prism  (without  the  entrance  channel  -  1.4  million  cubic 
metres) 

The  tidal  prism  was  calculated  by  multiplying  the  surface  area  of  the  reser¬ 
voir  without  the  entrance  channel  (2.5  million  square  metres)  by  the  average 
tidal  range. 

11.  Based  on  these  rough  estimates  of  replacement  times,  it  appears 
that  the  flushing  action  of  withdrawals  dominates  in  the  main  stem  of  the 
reservoir  between  the  entrance  and  the  withdrawal  points.  Dye  tests  reported 
in  this  Part  and  numerical  simulations  reported  in  Part  III  evaluated  flushing 
in  this  area.  In  other  back,  off-channel  areas  of  the  reservoir,  especially 
Foster  Creek  and  the  upper  part  of  Back  River,  the  largest  flushing  action 
comes  from  water  level  fluctuations,  also  evaluated  in  Part  III. 

Dve  Test  Procedures 

12.  Three  dye  tests  were  performed  by  USGS  starting  29  June,  9  August, 
and  30  August  1989  and  carried  out  over  3-day  periods.  The  first  test  was 
used  to  establish  procedures,  but  the  quantity  of  dye  injected  was  insuffi¬ 
cient  to  trace  with  confidence.  The  second  and  third  tests  used  34  and  42  kg 
of  Rhodamine-WT  fluorescent  dye,  respectively.  The  third  test  produced  the 
greatest  quantity  of  usable  dye  concentration  data. 

13.  Special  monitoring  sites  were  operated  specifically  for  this  study 
by  USGS  within  Bushy  Park  Reservoir  and  are  shown  in  Figure  1,  Insert  A.  Not 
shown  are  Ites  1-3,  located  on  the  Cooper  River.  Water  level  gage  reading 
and  gage  datums  from  these  sites  were  provided  to  WES  by  USGS.  The  datum 
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provided  for  site  8  appeared  to  be  in  error,  and  was  adjusted  by  0.3  m  based 
on  comparisons  to  nearby  gages. 

14.  In  addition  to  water  level  measurements,  USGS  measured  discharges 
half-hourly  at  two  locations  (near  sites  5  and  8)  along  the  existing  entrance 
canal  over  complete  tidal  cycles.  Supplemental  data  included  hourly  and  daily 
average  inflows  at  Pinopolis  provided  by  the  Charleston  District,  and  SCE&G 
provided  daily  average  withdrawals. 

15.  Flows  were  measured  by  integrating  current  speeds  over  the  two 
cross  sections  in  the  entrance  canal.  Cross-sect? on  depths  were  surveyed,  and 
a  tag  line  cable  tied  across  at  each  range.  Two  boats  were  used  at  each  range 
to  measure  currents  half-hourly  at  0.2  and  0.8  of  the  instantaneous  depth  and 
3— m  intervals  along  the  tag  line.  Gurley  current  meters  with  audio  output 
were  used  by  USGS  to  measure  current  speeds.  Raw  flow  data  were  processed  by 
USGS,  and  results  for  flow  versus  time  for  the  two  ranges  were  provided  to 
WES. 

16 .  Flow  measurements  began  near  uhe  end  of  ebb  tide  cycles .  At  the 
first  indication  that  the  flow  had  completely  turned  to  flood,  dye  injection 
began.  Dye  was  injected  at  site  4  over  5-  to  10-min  periods  from  a  boat  pass¬ 
ing  slowly  back  and  forth  across  the  canal.  A  long  bright  red  cloud  was  thus 
immediately  formed  in  the  canal  and  rapidly  expanded  into  the  reservoir. 

17.  Three  boats  were  used  for  hourly  daylight  dye  sampling.  After  the 
first  day  of  each  test,  two  boats  were  used  to  sample  dye.  Samples  were  taken 
at  surface,  middepth,  and  bottom.  Surface  and  bottom  sampling  points  were 
0.61  m  below  the  surface  and  0.61  m  up  from  the  bottom.  Small  water  bottles 
were  used  to  collect  samples.  At  night,  automatic  samplers  were  used  to 
collect  samples  at  sites  6,  10,  and  15. 

18.  Samples  were  analyzed  for  dye  concentration  in  the  field  using  a 
Turner  Designs  Model  10  fluorometer.  Dye  sample  times,  locations,  and  concen¬ 
trations  were  provided  to  WES . 

Results  and  Discussion 

19.  Selected  data  from  the  dye  test  starting  30  August  1989  are  pre¬ 
sented  as  follows.  Figure  2  shows  water  level  data  for  the  entrance  canal  and 
the  reservoir  over  the  3-day  period  of  the  third  dye  test.  Stages  in  Figure  2 
and  succeeding  figures  are  referred  to  the  National  Geodetic  Vertical  Datum 
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(NGVD) .  Site  6  had  an  hourly  recording  interval,  while  other  sites  had  half- 
hourly  intervals.  Figure  3a  shows  the  data  from  sites  8  and  15  plotted  along 
with  the  average  stage  for  the  reservoir.  Figure  3b,  the  stage  difference 
between  site  8  and  site  15,  shows  that  differences  can  be  over  0.15  m  during 
flood  tide  phases,  and  over  0.1  m  during  ebb  tidal  phases.  These  data  were 
interpolated  to  10-min  intervals  using  a  spline.  Stage  differences  should  be 
regarded  as  being  relative  since  datums  are  unverified. 

20.  Raw  dye  data  consist  of  concentrations  at  various  points  and  times. 
Data  were  divided  into  two  groups:  sites  3-8  for  the  entrance  canal  and 
sites  9-15  for  the  reservoir  proper.  Figure  4  shows  surface,  middepth,  and 
bottom  concentrations  measured  30  August-1  September  1939  for  sites  9-15. 

Mi  \depth  samples  were  not  always  collected.  Concentrations  are  reported  in 
dye  units  based  on  fluorometer  response  to  calibration  dilutions. 

21.  No  dye-decay  or  background  corrections  were  applied  to  the  data. 

Dye  photodecomposition,  chemical  degradation,  and  loss  by  adsorption  to  sedi¬ 
ments  were  deemed  too  small  to  be  important  for  these  dye  tests.  However, 
there  might  have  been  dye  loss  by  aquatic  plant  absorption,  and  uptake  tests 
were  performed  to  gage  this  effect. 

22.  Tests  on  the  uptake  of  Rhodamine-WT  by  water  hyacinths,  similar  to 
those  which  surround  the  reservoir  to  depths  of  about  9  ft,  were  conducted  at 
WES.  Four  small,  well-rooted  plants  were  crowded  into  three  tanks,  and  three 
additional  tanks  were  used  for  controls.  Dye  concentrations  in  the  hyacinth 
tanks  were  reduced  by  9  percent  during  the  first  day,  and  by  12  percent  during 
the  first  '  day?  ^mpared  to  control  tanks.  Plant  uptake  of  dye  is  thought  to 
depend  strongly  on  the  growth  state  of  the  plants,  which  is  difficult  to 
assess  both  for  the  laboratory  and  the  field.  Hyacinths  dominate  in  shallow 
areas  of  the  Bushy  Park  Reservoir,  while  deep  areas  of  greater  volume  and  flow 
are  free  from  the  plants.  Because  of  the  difficulties  applying  the  uptake 
test  results  to  the  field,  no  direct  corrections  were  made  to  the  dye  data; 
but  the  uptake  rates  were  used  to  establish  a  qualitative  sense  of  results 
accuracy . 

23.  Canal  flow  data  for  the  second  and  third  tests  are  presented  in 
Part  III.  Daily  freshwater  releases  from  the  Pinopolis  Dam  near  the  time  of 
the  third  dye  test  were  as  follows: 


Date 

Mean  Discharge 
cu  m/sec 

8/26 

73.6 

8/27 

110.9 

8/28 

251.7 

8/29 

202.5 

8/30 

88.9 

8/31 

69.6 

9/1 

93.1 

The  exact  total  withdrawal  rate  is  not  known  for  the  period  of  the  dye  tests, 
but  data  obtained  through  the  Bushy  Park  Vater  User's  Association  indicated 
that  the  largest  withdrawal  (24.1  cu  m/sec  daily  average)  was  at  capacity 
during  the  third  test. 


Analysis  of  Dve  Data 


Tidal  prism  model  description 

24.  A  simple  model  was  used  to  analyze  the  dye  data  from  the  third  dye 
test  and  to  evaluate  reservoir  flushing.  The  model  calculates  the  mass  and 
concentration  of  a  tracer  over  tidal  cycles  based  on  conservation  principles. 
The  model  is  similar  to  a  previous  tidal  basin  flushing  model  that  has  shown 
good  agreement  with  physical  hydraulic  model  results  (Callaway  1981) .  With¬ 
drawal  losses  were  added  to  the  model. 

25.  The  assumption  is  made  that  the  reservoir  is  well  mixed.  The  tidal 
prism  Vp  -  Vh  -  V,  where  V  is  volume,  and  the  subscripts  h  and  i  refer 
to  tidal  high  and  low  water  stages,  respectively.  Then  Vp  -  A„R  where  Aq 

is  the  reservoir  surface  area  and  R  is  the  tidal  range;  and  (Vp  -  withdraw¬ 
al)  is  the  ebb  flow  volume  through  the  entrance  canal.  If  the  reservoir  with¬ 
drawal  rate  is  Q„  and  the  tidal  period  is  T  ,  then  the  mass  of  tracer  M 
in  the  reservoir  at  the  end  of  ebb  is : 

M,  -  \  -  C*  R  ♦  Q„  *J  (3) 

where  Mi,  and  Cj,  are  the  mass  and  concentration,  respectively,  at  the 
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previous  high  water.  The  concentration  of  tracer  does  not  change  by  dilution 
during  the  ebb  tide  (C,  -  previous  Ch) .  During  the  flood  tide,  the  mass  of 
tracer  in  the  reservoir  is  not  changed  by  the  tidal  flow,  but  is  reduced  by 
the  withdrawal.  Thus,  at  the  end  of  the  flood  tide: 


Mh  -  M,  -  C,  Q*  * 


(4) 


and 


Ch  -  Mh  /  Vh 


(5) 


The  model  is  used  to  evaluate  concentrations  and  masses  of  tracer  over  a 
series  of  tidal  high  and  low  waters  in  steps. 

Dye  data  averaging 

26.  Data  plots  such  as  Figures  4a-c  showed  that  the  dye  was  not  uni¬ 
formly  distributed  in  the  reservoir,  and  that  concentrations  had  considerable 
scatter  in  time  as  well.  For  comparison  to  model  values,  dye  data  were  aver¬ 
aged  over  flood  and  ebb  tidal  phases.  Sampling  sites  were  located  along  the 
channel  axis  of  the  reservoir.  Data  were  therefore  averaged  in  space  by  first 
depth  averaging  and  then  weighing  according  to  the  site  depth.  Average  con¬ 
centrations  were  then  normalized  by  the  initial  observed  concentration  C0  . 
Rfeguisa,  .ftnalygeg 

27 .  Model  analyses  for  reservoir  flushing  were  performed  both  without 
withdrawals  and  with  an  assumed  28.3-cu  m/sec  withdrawal  rate.  Results  of  the 
two  model  analyses  and  for  the  averaged  dye  data  are  shown  in  Figure  5.  Dye 
test  results  should  be  compared  to  model  predictions  that  included  with¬ 
drawals.  The  flushing  shown  by  the  dye  tests  was  somewhat  greater  than  the 
model  predicted,  but  given  the  difficulties  inherent  in  dye  tests  and  the 
simplified  assumptions  of  the  model,  the  agreement  was  acceptable.  The  ini¬ 
tial  dip  exhibited  in  the  dye  data  was  disregarded.  The  flushing  times  for 

50  percent  reduction  in  initial  concentrations  were  interpolated  as  1.8  and 
2.4  tidal  cycles  for  the  dye  data  and  the  model  that  included  withdrawals, 
respectively.  The  model  flushing  time  for  50  percent  reduction  in  initial 
concentrations  with  no  withdrawals  was  not  reached  in  the  four  tidal  cycles 
for  which  data  were  available.  Using  an  average  tide  to  drive  the  model  after 
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four  tidal  cycles,  50  percent  reduction  was  reached  in  five  tidal  cycles. 

Thus  the  model  indicated  that  the  withdrawals  reduced  flushing  times  by  about 
50  percent.  These  calculated  flushing  times  are  consistent  with  the  rough 
estimates  of  paragraph  10. 
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PART  III:  TIDAL  MODEL  FLUSHING  TESTS 

Description  of  Model  Tests 

28.  The  two-dimensional  laterally  averaged  numerical  model  called  Fine¬ 
grained  Bed  Sediments  (FIBS)  was  used  for  this  study.  The  model  was  pre¬ 
viously  applied  to  predict  the  effect  of  the  1.5-m  Charleston  Harbor  channel 
deepening  on  shoaling  (Teeter  and  Pankow  1989).  However,  model  features  that 
stress  fine-grained  bed  sediments  were  not  employed  in  this  application. 

29.  The  model  numerically  solved  laterally  averaged  dynamic  differen¬ 

tial  equations  for  flow  continuity  and  horizontal  momentum  over  the  interval 
between  the  water  surface  and  the  channel  bottom  and  along  the  length  of  the 
channel.  The  horizontal  momentum  equation  included  nonlinear  advection  and 
quadratic  bottom  friction  terms.  FIBS  solves  equations  over  finite  elements 
using  a  method  of  weighted  residuals .  The  flow  domain  was  discretized  as  a 
series  of  elements  with  15  computational  points,  configured  as  five  roughly 
horizontal  rows  of  three  nodes.  Each  node  consisted  of  a  vertical  column  of 
five  computational  points.  Breadth  at  each  node  varied  vertically  as  a  quad¬ 
ratic  function  of  depth.  Model  equations  were  solved  explicitly  over  elements 
and  with  time.  Appendix  A  presents  a  brief  description  of  the  model. 

30.  The  numerical  mesh  covered  Charleston  Harbor  and  the  Cooper  River 

from  the  jetties  (river  mile  0)  to  the  Pinopolis  Dam  (river  mile  57)  at 
1.6-km  node  or  cross-section  spacing.  The  mesh  included  12.9  km  of  the  Wando 
River  and  25.7  km  of  the  East  Branch  Cooper  River,  also  at  1.6-km  node  spac¬ 
ing.  Bushy  Park  Reservoir  was  represented  in  the  model  at  0.8-km  node  spacing 

from  the  entrance  to  the  dam  and  also  4  km  of  Foster  Creek. 

31.  The  model  was  operated  by  specifying  water-surface  elevations  at 
the  ocean  (river  mile  0)  and  inflows  at  Pinopolis.  Optionally,  withdrawals 
were  made  in  Bushy  Park  Reservoir  at  two  points,  near  the  SCE&G  plant  (the 
largest  withdrawal)  8.8  km  from  the  entrance  and  near  the  dam,  as  described 
later . 

Model  Operation 


fmotm  &&A 

32.  Water  levels  are  monitored  along  the  length  of  the  Cooper  River  by 
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USGS.  Figure  1  shows  the  location  of  stations  used  for  model  verification. 

In  addition,  special  monitoring  station  sites  were  operated  specifically  for 
this  study  by  USGS.  Those  sites  are  located  within  Bushy  Park  Reservoir,  as 
shown  in  Figure  1,  Insert  A.  In  addition  to  water  level  measurements,  USGS 
measured  discharges  half-hourly  at  two  locations  (sites  5  and  8)  along  the 
existing  entrance  canal  over  two  separate  tidal  cycles.  Because  vertical 
velocity  measurements  are  difficult  to  obtain,  no  vertical  velocity  data  were 
available  with  which  to  verify  the  model. 

33.  Hourly  and  daily  average  inflows  at  Pinopolis  Dam  were  provided  by 
the  Charleston  District.  SCE&G  provided  daily  average  withdrawals. 
Verification  procedure 

34.  Ocean  boundary  specifications  were  developed  from  hourly  Customs 
House  water  level  data,  fit  with  a  cubic  spline.  Water  level  values  were 
extracted  at  172-sec  intervals  (the  model  time-step),  and  multiplied  by  0.9. 
The  resulting  boundary  condition  was  found  to  produce  reasonable  water  level 
results  at  the  Customs  House  (river  mile  9.5)  in  the  model.  Freshwater  inflow 
specifications  were  developed  from  hourly  Pinopolis  Dam  discharge  readings 
smoothed  by  locally  weighted  regression. 

35.  The  model  was  verified  to  two  data  periods.  The  first  data  period, 
11-14  June  1988,  corresponded  to  the  onset  of  a  salinity  intrusion  event  exam¬ 
ined  in  Part  IV.  This  data  period  was  used  to  verify  the  Cooper  River  water 
levels  from  the  Customs  House  to  Pimlico. 

36.  For  the  first  data  period,  the  model  was  first  spun  up  from  a 
quiescent,  uniform  water  level  to  a  dynamic  condition  for  two  tidal  cycles. 
Results  from  those  calculations  were  used  to  restart  the  model  for  a  second 
set  of  calculations  lasting  five  tidal  cycles.  The  final  three  cycles  of  this 
run  were  used  for  comparison  to  field  data. 

37.  The  second  data  period,  29-31  August  1989,  corresponded  to  an 
intensive  data  collection  period  by  USGS  as  described  in  Part  II.  The  data 
considered  during  this  period  came  from  Bushy  Park  Reservoir  entrance  to  the 
dam.  The  second  data  period  also  included  flow  measurements  within  the 
entrance  channel  at  sites  5  and  8.  The  procedure  for  operating  the  model  was 
similar  to  that  used  for  the  first  data  period  described  earlier. 

Elen  tests 

38.  The  proposed  relocation  route  was  used  for  all  plan  tests  (Fig¬ 
ure  1,  Insert  A).  The  first  plan  tested  used  the  average  cross-section 
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dimensions  of  the  existing  canal.  Subsequent  plan  tests  had  cross-section 
dimensions  covering  a  range  to  reproduce  the  existing  tidal  ranges  and 
entrance  canal  flows.  Plan  tests  used  tidal  and  inflow  conditions  from  the 
second  data  period  (29-31  August  1989). 

Canal  Test  Results  and  Discussion 


Model  verification 

39.  Figure  6  shows  comparisons  between  model  and  prototype  water  levels 
for  the  first  data  period  starting  12  June  1988 .  Note  that  data  shown  in 
Figure  6  for  Mempkin  were  actually  from  Pimlico,  2  miles  downstream.  Mempkin 
tide  gaging  was  only  recently  initiated  and  does  not  extend  to  this  data 
period.  Comparisons  of  Mempkin  and  Pimlico  water  levels  for  1-7  May  1989 
showed  water  level  fluctuations  to  be  almost  identical,  with  the  root  mean 
squares  (RMS)  about  the  means  differing  by  only  0.003  m.  Figure  7  shows  part 
of  that  comparison.  It  is  difficult  to  judge  phases  using  hourly  data,  but  it 
appears  that  the  two  stations  are  nearly  in  phase.  Because  of  the  importance 
of  Mempkin  as  the  interception  point  for  the  proposed  canal,  model  comparisons 
were  made  at  this  point  using  data  from  the  Pimlico  gage. 

40.  It  was  difficult  to  establish  accurate  descriptions  for  shallow, 
off-channel  areas,  especially  at  and  above  the  Tee  on  both  branches  of  the 
Cooper  River.  Many  of  these  areas  were  formerly  rice  ponds  and  partially 
leveed  from  the  river.  These  areas  act  as  off-channel  storage  for  flow  but  do 
not  contribute  directly  to  momentum  of  the  flow.  Satellite  images  and  air 
photos  were  useful  in  defining  top  width  near  and  below  the  Tee,  but  areas 
above  Mempkin  were  not  well  known. 

41.  Pinopolis  inflow  is  usually  zero  for  about  8  hr  a  day,  and  other¬ 
wise  varies  widely.  The  model  responded  more  rapidly  than  the  prototype  to 
rapid  changes  in  Pinopolis  inflow.  This  was  the  reason  for  smoothing  hourly 
inflow  data.  The  effect  can  be  seen  in  Figure  6  at  day  14.2  when  low  water 
was  exaggerated  in  the  model  by  falling  inflows  from  Pinopolis.  Previously 
mentioned  uncertainties  in  off-channel  storage  were  undoubtedly  responsible 
for  this  effect. 

42.  Also  note  that  DuPont  Chemical  water  levels  (Figure  6)  were  influ¬ 
enced  by  conditions  within  the  reservoir,  as  will  be  discussed  later.  Over¬ 
all,  verification  to  the  first  data  period  was  considered  good. 
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43.  Verification  to  the  second  data  period  starting  29  August  1989 
encountered  uncertainties  in  specifying  reservoir  withdrawals.  First  efforts 
used  constant  withdrawals  based  on  withdrawal  capacity,  but  were  not  success¬ 
ful  at  reproducing  water  levels.  Further  checks  with  the  Bushy  Park  Water 
Users  Association  and  with  SCE&G  indicated  that  withdrawals  are  variable  dur¬ 
ing  an  average  day.  The  daily  average  withdrawal  from  SCE&G  during  this 
period  was  24  cu  m/sec . 

44.  USGS  flow  data  for  the  second  data  period  complicated  the  with¬ 
drawal  question  by  showing  a  mean  or  net  flow  of  -0.79  cu  m/s ec,  with  the 
negative  sign  indicating  the  seaward  direction  (Figure  8).  The  flow  data  were 
not  in  agreement  with  the  reported  daily  withdrawal,  indicating  that  the  with¬ 
drawal  measured  during  the  tidal  cycle  was  very  small  or  that  it  was  offset  by 
a  decrease  in  Cooper  River  water  level  over  that  tidal  cycle.  During  another 
intensive  survey  performed  on  9  August  1989,  the  tidal  flow  record  did  exhibit 
the  expected  offset  of  about  31  cu  m/sec  corresponding  to  the  withdrawal 
amount  (Figure  9) . 

45.  Hourly  withdrawals  were  not  known.  To  examine  the  effects  of 
withdrawals  on  water  levels  and  flows,  three  withdrawal  scenarios  were  tested: 
no  withdrawal,  capacity  withdrawals,  and  a  variable  withdrawal  averaging  the 
capacity.  The  prototype  data  generally  fell  between  the  constant  no¬ 
withdrawal  and  capacity-withdrawal  cases  (Figures  10  and  11).  Variable  with¬ 
drawals  were  found  to  either  increase  or  decrease  water  level  fluctuations 
depending  on  timing.  Several  variable  withdrawal  schedules  were  tested. 

Figure  12  shows  the  scenario  in  which  withdrawals  were  reduced  to  8  cu  m/sec 
between  2200  hr  and  0600  hr,  and  the  withdrawals  during  the  rest  of  the  day 
were  such  that  the  daily  average  withdrawal  was  31  cu  m/sec.  This  withdrawal 
scenario  was  used  for  all  subsequent  plan  testing. 

46.  While  the  effects  of  withdrawals  can  be  seen  in  the  mean  water 
levels,  model  tests  of  these  withdrawal  scenarios  showed  that  the  average  tide 
range,  as  calculated  from  the  root  mean  square  deviations  of  water  levels  from 
the  mean,  was  insensitive  to  withdrawal  variations.  For  example,  the  average 
tide  ranges  for  site  15  shown  in  Figures  10-12  are  within  0.003  m  for  the  last 
three  full  tidal  cycles.  Thus,  even  though  water  levels  were  not  reproduced 
in  the  model  exactly,  the  model  verification  to  the  second  data  period  was 
considered  to  be  good. 
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Plan  tests 

47.  The  first  plan  of  the  lengthened  canal  tested  used  the  existing 
canal  cross  sections.  The  results  are  displayed  in  Figure  13.  Model  results 
critical  to  the  tidal  flushing  of  the  entire  reservoir  and  surrounding  tidal 
marshes  were  the  water  level  fluctuations  at  sites  8  and  15  (the  two  ends  of 
the  reservoir),  and  the  flows  in  the  entrance  canal.  Model  results  from  the 
last  three  complete  tidal  cycles  were  used  to  compute  water  level  and  flow 
variances,  and  these  were  then  converted  to  average  tidal  ranges  and  flow 
variabilities.  Results  of  model  base  and  plan  tests  are  shown  in  the 
following  tabulation: 


Three-Tidal-Cvcle  Average  Values 


Canal 

Width 

m 

Condition 

Depth 

m 

Tide  Range 

Site  8 

m 

Tide  Range 
Site  15 

m 

Flow  Variability 
RMS,  Site  4 
cu  m/sec 

Existing  Canal 

55 

4.2 

0.445 

0.434 

81.8 

Lengthened 

Canal 

55 

4.2 

0.366 

0.379 

63.6 

65 

5 

0.380 

0.391 

66.8 

75 

5 

0.394 

0.408 

67.5 

95 

5 

0.404 

0.418 

70.9 

65 

6.5* 

0.465 

0.483 

88.2 

65 

6** 

0.442 

0.461 

83.0 

*  Includes  deepened  old  canal  +  2  m  and  channel  at  site  8  +  1  m. 

**  Includes  deepened  old  canal  +  1.5  m  and  channel  at  site  8  +  1.5  m. 

The  effect  of  canal  lengthening  (with  existing  cross  section)  was  greater  on 
flows  (22  percent  reduction)  than  on  water  level  fluctuations  (16  percent 
reduction) .  Even  though  the  95-m-wide  by  5-m-deep  planned  canal  had  a  much 
larger  cross  section  (293  sq  m)  than  the  existing  canal  (164  sq  m) ,  tidal 
ranges  and  flows  were  not  improved  much.  (It  was  assumed  that  it  was  feasible 
only  to  deepen  the  existing  canal,  since  widening  would  involve  bridge  and 
levee  modifications.) 

48.  Figure  14  shows  model  results  for  the  most  successful  canal  cross- 
section  design:  65  m  wide  by  6  m  deep.  Figure  15  shows  typical  cross- 
sections  for  the  proposed  new  canal,  remaining  portions  of  the  existing 
entrance  canal,  and  the  revised  reservoir. 


Withdrawal  Scheme  Test 

49.  A  withdrawal  scheme  to  augment  flushing  with  the  smaller,  existing 
entrance  canal  cross  section  was  tested.  Results  of  model  tests  using  the 
existing  canal  dimensions  for  the  proposed  canal  cross  section  indicated  that 
the  tidal  component  of  flushing  decreased  by  22  percent  compared  to  the  exist¬ 
ing  condition  (paragraph  47).  Using  half  the  existing  flow  variability 

(40.9  cu  m/sec)  as  the  average  ebb  flow  over  a  tidal  cycle,  the  tidal  flushing 
component  can  be  estimated  to  decrease  by  about  9  cu  m/sec  for  the  case  of  the 
lengthened  canal  using  the  existing  cross  section.  Thus,  a  supplemental  with¬ 
drawal  of  9  cu  m/sec  would  compensate  for  the  difference  in  the  cross-section 
area  between  the  proposed  and  existing  canal  dimensions,  neglecting  any 
effects  of  the  supplemental  withdrawal  on  tides . 

50.  To  test  this  result,  the  model  was  rerun  for  the  following  condi¬ 
tions:  base  (existing),  lengthened  canal  with  existing  cross  section,  and 
lengthened  canal  with  existing  cross  section  including  a  9-cu  m/sec  supple¬ 
mental  withdrawal  at  the  Back  River  Dam.  The  model  was  spun  up  with 
127-cu  m/sec  inflow  at  Pinopolis  and  31-cu  m/sec  Bushy  Park  withdrawal  for 
three  tidal  cycles.  Model  results,  plotted  as  cumulative  flushing  volume 
(withdrawal  and  ebb  tidal  flow)  in  Figure  16  for  two  tidal  cycles,  confirm 
that  a  9-cu  m/sec  supplemental  withdrawal  compensates  for  the  reduced  tidal 
flushing.  The  supplemental  withdrawal  allowed  the  lengthened  canal  flushing 
volume  to  catch  up  with  the  base  condition  at  the  end  of  the  flood  tidal 
phase.  The  model  run  was  slightly  less  than  two  complete  tidal  cycles,  so 
that  this  occurred  only  once  (Figure  16). 

51.  Data  from  the  model  runs  were  also  used  to  explore  the  feasibility 
of  a  gravity  flow  outlet  by  comparing  the  water-surface  elevations  across  the 
Back  River  Dam.  Figure  17  shows  two  tidal  cycles  of  water  levels  on  the  river 
side  and  on  the  reservoir  side,  including  the  same  three  reservoir  conditions 
described  in  paragraph  50.  Tidal  phase  and  amplitude  were  altered  at  the 
Bushy  Park  side  of  the  dam  by  the  lengthened  canal,  and  were  almost 

180  degrees  out  of  phase  across  the  dam.  The  head  difference  across  the  dam 
was  as  much  as  a  metre  for  the  supplemental  withdrawal  test,  as  shown  in  Fig¬ 
ure  18.  The  supplemental  withdrawal  was  used  to  demonstrate  (rather  crudely) 
the  effect  of  reservoir  drawdown.  Rough  calculations  indicate  that  an  effi¬ 
cient  weir  of  20-  to  30-m  width  (with  flap  gate)  might  be  sufficient  to 


19 


produce  the  desired  9-cu  m/sec  average  supplemental  withdrawal  from  the 
reservoir. 

52.  Preliminary  design  of  the  proposed  entrance  canal  using  recent 
Charleston  District  topographic  information  and  cross  sections  at  61-m  inter¬ 
vals  indicated  the  oversized  canal  would  require  excavation  of  2.08  million 
cubic  metres  as  compared  to  1.42  million  cubic  metres  for  the  existing  canal 
dimensions.  Therefore,  substituting  a  pumped  or  weir  withdrawal  for  the  over¬ 
sized  canal  would  save  0.65  million  cubic  yards  of  excavation. 

53.  The  substitution  of  withdrawal  for  tidal  volumes  to  achieve  flush¬ 
ing  would  be  most  effective  along  the  main  stem  of  the  reservoir.  The  rela¬ 
tive  effects  of  the  tidal  and  withdrawal  components  may  not  be  equivalent  in 
backwater  areas  of  Foster  Creek  and  Back  River,  as  noted  earlier. 


20 


PART  IV:  SALINITY  CONDITIONS 


54.  Data  from  intrusion  events  were  used  to  drive  an  analytic  model  of 
salinity  intrusion.  Under  reduced  freshwater  inflow  conditions,  ocean  chlo¬ 
rides  intrude  farther  upstream  in  the  Cooper  River.  The  reservoir  has  experi¬ 
enced  occasional  ocean  chloride  intrusion  events  since  rediversion,  and  values 
have  reached  as  high  as  100  ppm.  It  is  these  intrusion  events  that  have  led 
the  Charleston  District  to  continue  to  consider  structural  measures  to  protect 
Bushy  Park  water  quality  and  inflow  scheduling  options  to  improve  flexibility 
of  the  hydroelectric  plant  operated  by  the  State  of  South  Carolina. 

55.  A  Charleston  District  alert  system,  operated  by  USGS,  senses 

chloride  and  water  levels.  Based  on  releases  from  the  Pinopolis  Dam  over  the 

previous  several  days,  the  system  attempts  to  anticipate  conditions  conducive 

to  intrusion  events.  The  system  automatically  orders  supplemental  freshwater 
releases  from  the  Pinopolis  Dam  when  necessary  to  protect  Bushy  Park  water 
quality.  However,  because  of  the  time  lag  required  to  halt  a  chloride  intru¬ 
sion  event,  supplemental  releases  are  only  partially  effective.  Water  level 

and  specific  conductance  fluctuations  along  the  length  of  the  Cooper  River  are 
shown  in  Figure  19  for  a  typical  intrusion  event. 

Canal  Relocation  Effectiveness 

56.  The  effectiveness  of  the  canal  relocation  was  gaged  by  calculating 
salinities  for  previous  intrusion  events.  Field  data  and  an  analytic  model 
were  used.  The  analytic  model  is  a  one-dimensional  closed-form  solution  of  an 
estuarine  convection-diffusion  equation  which  includes  tidal  convection 
effects  (tide  range  and  tidal  velocity  at  the  estuary  mouth),  freshwater  in¬ 
flow,  and  water  depth  (Headquarters,  US  Army  Corps  of  Engineers,  1991).  The 
model,  which  assumes  uniform  distributions  of  salinity  and  velocity  and  con¬ 
stant  depth  and  width  over  cross  sections,  calculates  longitudinal  salinity 
profiles  for  given  sets  of  conditions.  The  model  was  used  to  calculate  salin¬ 
ities  at  the  proposed  entrance  location  for  intrusion  events.  These  values 
were  corrected  for  the  hydraulic  effect  of  the  proposed  canal  and  converted  to 
chlorides . 

57.  Existing  data  were  used  to  project  ocean  chloride  intrusion  at  the 
proposed  entrance  canal  location.  The  critical  stations  are  DuPont  Chemical 
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and  Mempkin  for  the  existing  and  proposed  canal  sites,  respectively.  Seven 
extreme  intrusion  events  were  identified  and  used  for  analyses.  Selected  data 
were  used  as  model  input  to  calculate  the  longitudinal  distributions  of 
salinity,  and  other  data  were  used  for  model/prototype  comparison. 

58.  Conversions  between  specific  conductance  (field  data),  salinity 
(model  data) ,  and  chloride  concentration  (water  quality  criteria)  were 
required.  Specific  conductance  values  from  the  USGS  gages  were  converted  to 
salinity  using  standard  oceanographic  methods.  These  methods  are  inaccurate 
at  low  salinities,  where  ratios  between  ionic  species  vary,  so  another  method 
was  required.  Chloride  concentration  is  customarily  used  to  assess  water 
quality  conditions  in  Bushy  Park  Reservoir.  Specific  conductance  SC  and 
chloride  Cl  values  (134  data  points)  from  previous  studies  in  the  Cooper 
River  have  been  used  to  develop  an  empirical  regression  relationship: 


Cl  =  0.2398  *  SC  -  4.2636 


(6) 


where  Cl  is  in  ppm  and  SC  is  in  microsiemens  (micro-S)  per  cm.  This  rela¬ 
tionship  is  applicable  to  river  and  estuarine  waters  with  SC  less  than  about 
2,000  micro-S  per  cm,  and  has  a  standard  error  of  ±32  ppm  Cl  associated  with 
it.  Chlorinities  C  in  ppm  can  be  calculated  from  salinities  S  in  ppm 
using  the  oceanographic  expression: 


C  =  S/1.80665 


(7) 


It  was  found  that  values  of  C  and  Cl  were  closely  related,  and  that 


Cl  *  0.736  *  C 


(8) 


could  be  used,  along  with  the  other  expressions  described,  to  compare  model 
and  field  results. 

59.  The  model  requires  salinity  values  at  a  minimum  of  three  locations 
as  input.  The  ocean  source  salinity  was  specified  as  34  ppt  in  the  model  at 
the  mouth  of  the  estuary.  Salinity  values  from  Mobay  Chemical  and  Pimlico 
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stations  were  used  for  model  input  at  the  other  two  locations.  Values  were 
calculated  for  other  stations,  including  the  proposed  relocation  site  for  the 
entrance  canal  at  Mempkin. 

60.  The  model  calculations  of  maximum  salinity  at  the  relocation  site 
(Mempkin)  were  corrected  to  account  for  alteration  in  tidal  excursions  in  the 
vicinity  of  the  entrance  canal.  The  flow  in  the  Cooper  River  between  the 
existing  and  proposed  entrance  canals  would  be  altered  by  the  canal  reloca¬ 
tion.  During  flood  tide,  flow  that  had  previously  entered  Bushy  Park  would 
continue  upstream  and  feed  into  the  relocated  canal.  The  correction  was  cal¬ 
culated  as  the  product  of  the  increased  tidal  excursion  length  and  the  longi¬ 
tudinal  salinity  gradient  in  the  Cooper  River.  The  excursion  length  increase 
was  estimated  as  2.78  km  for  the  tide  ranges  and  tidal  prisms  expected  to 
occur  during  intrusion  events. 


Information  Sources 

61.  USGS  maintains  for  the  Charleston  District  a  number  of  r^al-time 
conductivity  and  water  level  stations  (including  other  water  quality  parame¬ 
ters)  along  the  Cooper  River.  Special  short-term  stations  have  also  been 
installed  and  operated  in  certain  areas  of  interest.  Stations  generally 
record  data  internally  at  15-min  intervals  and  ^  _ .ad  data  every  4  hr  via 
satellite  link  to  the  USGS  computer  at  Columoia,  SC.  Stations  used  in  these 
study  tasks,  located  in  Figure  1,  include  Pimlico,  DuPont  Chemical,  Goose 
Creek,  Mobay  Chemical,  Customs  House,  and  Army  Depot.  The  DuPont  Chemical 
station  is  located  along  the  existing  Bushy  Park  entrance  canal.  Figure  1 
also  shows  Mempkin,  the  site  of  the  proposed  canal  relocation.  Conductivity 
and  water  level  data  were  inspected  and  selectively  downloaded  from  the  USGS 
computer.  Downloaded  sensor  data  were  averaged  by  hour  and  separated  by 
category. 


Results  and  Discussion 

62.  Table  1  lists  daily  maximum  instantaneous  conductance  values  at 
DuPont  Chemical  for  days  bracketing  the  seven  worst  intrusion  events.  The 
seven  worst  days,  listed  in  the  following  tabulation,  were  used  in  the  model 
analysis.  Figure  20  shows  three  example  longitudinal  salinity  distributions 
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Daily  Maximum  Chlorides,  ppm 


Event 

Date 

Data 

Source 

Mob ay* 
Chemical 

Goose 

DuPont 

Chemical 

Pimlico* 

MemDkin** 

6/24/86 

Field 

1,930 

229 

56 

40 

NA 

Model 

1,930 

279 

72 

40 

27 

7/21/86 

Field 

3,076 

242 

101 

55 

NA 

Model 

3,076 

437 

104 

55 

35 

12/22/87 

Field 

1,000 

101 

74 

31 

NA 

Model 

1,000 

168 

51 

31 

26 

6/04/88 

Field 

2,909 

217 

90 

37 

NA 

Model 

2,909 

355 

75 

37 

30 

6/25/88 

Field 

5,553 

315 

69 

36 

NA 

Model 

5,553 

563 

86 

36 

3 

7/28/88 

Field 

1,823 

185 

118 

37 

NA 

Model 

1,823 

262 

68 

37 

31 

7/29/88 

Field 

1,585 

537 

76 

37 

NA 

Model 

1,585 

241 

66 

27 

31 

*  Field  data  were  used  as  model  input. 

**  Corrected  for  tidal  excursion  length  (see  paragraph  60).  NA 
-  data  not  available. 


calculated  by  the  model.  Model  results  and  computed  Cl  values  for  field 
data  are  summarized  in  the  tabulation.  The  average  difference  between  field 
and  model-calculated  values  for  DuPont  was  9  ppm,  and  the  standard  deviation 
of  the  differences  was  24  ppm.  Data  from  water  years  1981-1983  previously 
established  the  background  chloride  concentration  in  the  Cooper  River  inflow 
to  vary  monthly  between  14  and  37  ppm,  and  average  26  ppm.  The  model  calcula¬ 
tions  show  that  conditions  at  the  proposed  entrance  canal  would  remain  near 
normal  (less  than  35  ppm  Cl)  during  similar  future  intrusion  events,  and  con¬ 
firm  that  the  canal  relocation  would  be  effective  in  eliminating  ocean 
chloride  intrusions  from  the  reservoir. 
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PART  V:  SUMMARY  AND  RECOMMENDATIONS 


63.  Dye  tests  were  performed  in  Bushy  Park  Reservoir  to  evaluate  reser¬ 
voir  flushing  under  the  present  reservoir  configuration  and  conditions. 
Selected  data  are  presented  and  analyzed  herein.  Flushing  was  found  to  be 
relatively  rapid  and  analysis  indicated  that  industrial  withdrawals  substan¬ 
tially  enhanced  this  process.  Observed  flushing  time  for  50  percent  replace¬ 
ment  was  1.8  tidal  cycles.  Results  were  representative  of  conditions  along 
the  main  axis  of  the  reservoir  where  dye  samples  were  collected  and 
withdrawals  are  made. 

64.  A  simple  tidal  flushing  model  with  reservoir  industrial  withdrawals 
incorporated  was  in  qualitative  agreement  with  field  dye  test  results.  With¬ 
drawals  substantially  improved  the  flushing  action  of  the  main  axis  of  the 
reservoir  between  the  entrance  and  the  withdrawal  points.  The  model  indicated 
about  a  50  percent  improvement  in  flushing  time  due  to  the  withdrawals.  Dye 
tests  carried  out  by  the  US  Geological  Survey  (USGS)  evaluated  the  rapid 
flushing  of  the  reservoir  main  axis  by  the  combined  action  of  tidal  exchanges 
and  withdrawals.  It  follows  that  in  other  back  areas  of  the  reservoir, 
especially  Foster  Creek  and  Back  River,  the  largest  flushing  action  comes  from 
water  level  fluctuations  which  were  the  focus  of  the  numerical  model  study 
performed. 

65.  The  effect  of  the  proposed  entrance  canal  on  tidal  flushing  of 
Bushy  Park  Reservoir  was  tested  using  a  numerical  model,  and  a  canal  cross 
section  that  would  maintain  present  flushing  was  developed  by  trial.  A  previ¬ 
ously  developed  two-dimensional  model  (Teeter  and  Pankow  1989)  computed  dy¬ 
namic  flow  velocities,  vertical  velocities,  and  water-surface  elevations  for 
two  verification  periods  and  for  various  entrance  canal  cross  sections.  The 
model  was  verified  to  tidal  volume  exchanges  and  water  levels  measured  by  USGS 
(vertical  velocity  calculations  were  not  verified) .  Observed  Charleston  Cus¬ 
toms  House  tides  and  Pinopolis  discharges  were  used  to  drive  the  model  bound¬ 
aries.  The  model  included  the  harbor  channel  configuration  (which  was  being 
deepened),  and  reservoir  withdrawals  as  they  existed  for  the  USGS  field  test. 
After  verification,  the  model  was  used  to  test  a  number  of  entrance  canal 
geometries  (width  and  depth)  for  the  proposed  relocation,  and  the  reservoir 
tidal  volume  exchanges  and  flushing  were  compared  to  existing  conditions. 

66.  Existing  tidal  flushing  can  be  maintained  in  Bushy  Park  Reservoir 
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by  constructing  an  upstream  entrance  canal  of  sufficient  size.  The  recom¬ 
mended  canal  cross  section  is  65  m  wide  by  6  m  deep,  with  side  slopes  of  about 
IV  on  3H,  and  a  cross-sectional  area  of  272  sq  m  at  mean  tidal  level.  In 
addition,  it  is  recommended  that  the  remaining  portion  of  the  existing  canal 
be  deepened  by  +1.5  m  to  5.65  m  over  a  length  of  3.2  km,  and  that  the  first 
0.8  km  of  the  reservoir  be  deepened  by  +1.5  m.  A  9-cu  m/sec  withdrawal  scheme 
was  tested,  and  results  indicated  that  such  a  pumped  or  gravity  weir  with¬ 
drawal  might  substitute  for  the  effect  of  the  canal  enlargement,  at  least 
along  the  main  stem  of  the  reservoir. 

67.  Prior  to  finalizing  the  entrance  canal  relocation  design,  it  is 
recommended  that  additional  cross  sections  be  taken  at  the  two  bridges 
crossing  the  existing  canal,  at  the  existing  entrance,  and  at  the  proposed 
entrance.  This  supplemental  information  will  be  used  to  verify  the  assump¬ 
tions  made  in  establishing  the  model  geometry.  It  is  also  recommended  that, 
prior  to  the  completion  of  the  new  canal  design,  an  analysis  should  be  made  of 
possible  sedimentation  within  the  proposed  canal.  In  addition,  previous  and/ 
or  new  soil  borings  should  be  evaluated  to  confirm  the  stability  of  the  IV  on 
3H  slope  used  in  this  study. 

68.  Recent  ocean  chloride  intrusion  data  were  used  to  confirm  that  the 
proposed  site  for  the  Bushy  Park  entrance  canal  would  effectively  reduce 
chloride  intrusions.  Analysis  of  previous  intrusion  events  indicated  that 
relocation  of  the  Bushy  Park  entrance  canal  to  the  proposed  site  near  Mempkin 
will  eliminate  the  salinity  intrusion  events  if  other  conditions  (weekly 
inflow  schedules,  for  instance)  remain  the  same. 


26 


REFERENCES 


Benson,  H.  A.  1976.  "Effects  of  40-ft  Charleston  Harbor  Project  on  Tides, 
Currents,  and  Salinities,"  Miscellaneous  Paper  HL-76-9,  US  Army  Engineer 
Waterways  Experiment  Station,  Vicksburg,  MS. 

_ .  1977.  "Charleston  Harbor  Navigation  Study,  South  Carolina,  Veri¬ 
fication  Tests,"  Miscellaneous  Paper  HL-77-4,  US  Army  Engineer  Waterways 
Experiment  Station,  Vicksburg,  MS. 

Benson,  H.  A.,  and  Boland,  R.  A.,  Jr.  1977.  "Dispersion  of  Proposed  Effluent 
Discharges  and  Saltwater  Intrusion  in  Cooper  River,"  Miscellaneous  Paper 
HL-77-14,  US  Army  Engineer  Waterways  Experiment  Station,  Vicksburg,  MS. 

Blumberg,  A.  F.  1977  (Mar).  "Numerical  Models  of  Estuarine  Circulation," 
Journal  of  the  Hydraulics  Division.  American  Society  of  Civil  Engineers. 

Vol  103,  No.  HY3 ,  pp  295-310. 

Blumberg,  A.  F. ,  and  Kantha,  L.  H.  1985.  "Open  Boundary  Condition  for  Circu¬ 
lation  Models , "  Journal  of  the  Hydraulics  Division.  American  Society  of  Civil 
Engineers .  Vol  111,  No.  2,  pp  237-255. 

Callaway,  R.  J.  1981  (May).  "Flushing  Study  of  South  Beach  Marina,  Oregon," 
Journal  of  Waterway.  Port.  Coastal,  and  Ocean  Engineering.  American  Society  of 
Civil  Engineers .  Vol  107,  No.  WW2,  pp  47-58. 

Finlayson,  B.  A.  1972.  The  Method  of  Weighted  Residuals  and  Variational 
Principles ,  Academic  Press,  New  York. 

Headquarters,  US  Army  Corps  of  Engineers.  1991  (25  Apr).  "Computation  of 
Salinity  Intrusion  by  One-Dimensional  Analysis,"  ETL  1110-8-7  (FR) ,  US 
Government  Printing  Office,  Washington,  DC. 

Munk,  W.  H. ,  and  Anderson,  E.  R.  1948.  "A  Note  on  the  Theory  of  the 
Thermocline,"  Journal  of  Marine  Research.  Vol  26,  pp  24-33. 

Sheng,  Y.  P.  1983.  "Mathematic  Modeling  of  Three-Dimensional  Coastal  Cur¬ 
rents  and  Sediment  Dispersion:  Model  Development  and  Application,"  Technical 
Report  CERC-83-2,  US  Army  Engineer  Waterways  Experiment  Station,  Vicksburg, 

MS . 

Smith,  T.  J.  1982.  "On  the  Representation  of  Reynolds  Stress  in  Estuaries 
and  Shallow  Coastal  Seas,"  Journal  of  Physical  Oceanography.  Vol  12,  pp  914- 
921. 

Smith,  T.  J.,  and  Takhar,  H.  S.  1979.  "On  the  Calculation  of  Width  Averaged 
Flow  Due  to  Long  Waves  in  an  Open  Channel,"  Journal  of  Hydraulic  Research. 

Vol  17,  pp  329-340. 

Teeter,  A.  M.  1989.  "Effects  of  Cooper  River  Rediversion  Flows  on  Shoaling 
Conditions  at  Charleston  Harbor,  Charleston,  South  Carolina,"  Technical  Report 
HL-89-3,  US  Army  Engineer  Waterways  Experiment  Station,  Vicksburg,  MS. 

Teeter,  A.  M. ,  and  Pankow,  W.  1989.  "Schematic  Numerical  Modeling  of  Harbor 
Deepening  Effects  on  Sedimentation,  Charleston,  South  Carolina,"  Miscellaneous 
Paper  HL-89-7,  US  Army  Engineer  Waterways  Experiment  Station,  Vicksburg,  MS. 

US  Army  Engineer  District,  Charleston.  1983  (March).  "Cooper  River  Rediver¬ 
sion  Project;  Special  Report:  Alternatives  for  Protection  of  Bushy  Park 
Reservoir,"  Draft  Report,  Charleston,  South  Carolina. 


27 


Wang,  D.  P.  1983.  "Two-Dimensional  Branching  Salt  Intrusion  Model,"  Journal 
of  Waterway.  Port.  Coastal,  and  Ocean  Engineering.  American  Society  of  Civil 
Engineers .  Vol  109,  No.  1,  pp  103-112. 

Wang,  D.  P.  1984.  "Mutual  Intrusion  of  a  Gravity  Current  and  Density  Front 
Formation,"  Journal  of  Physical  Oceanography.  Vol  14,  pp  1191-1199. 


28 


Table  1 


Year  Month  Days 


1986  June  18-26 


Conductivity 


1986  July  16-24 


1987  Dec.  19-25 


1988  May  30-31 


June  1-18 


PROPOSED 
CANAL  ROUTE 


PINOPOUS  DAM 


SITE  S 
SITES.  „ 


r'  C 

\  SITE  4  W. 

\  RICE  HOPE 


%  SITE  6  ( 

SITE  7* 


.os"'*  Qz~Pmjco  V 
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Figure  1.  Vicinity  map  and  sampling  site  locations 


levels  along  Bushy  Park  Reservoir  and  entrance  canal  starting 
30  August  1989  (the  third  dye  test) 
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Figure  3.  Water  levels  and  water  level  differences  for  Bushy  Park  Reservoir  starting 

30  August  1989 
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Figure  5.  Model  dye  washout  with  and  without  withdrawals,  and  averaged 
field  data  for  sites  9-15  starting  at  the  10: 00-16 :00-EST  ebb  tide  of 

30  August  1989 
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Figure  7.  Comparison  between  Mempkin  and  Pimlico  stages,  May  1989 


Figure  8.  Flows  at  sites  5  and  8  for  30  August  1989 
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results  for  second  data  period  starting  30  August  1989  with  no  reservoir 

withdrawals 
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a.  Entrance  canal  extension 


b.  Deepened  existing  canal 


Figure  15.  Cross  sections  for  the  proposed  lengthened 
entrance  canal  design  found  to  maintain  tidal  amplitude 
and  flushing  in  the  Bushy  Park  Reservoir 
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Figure  16.  Model  cumulative  flushing  volume  results  for  the  base, 
lengthened  canal,  and  lengthened  canal  with  9-cu  m/sec  supplemental 

withdrawal 


Figure  17.  Model  water  level  results  for  both  sides  of  the  Bushy  Park 
Dam  for  the  base,  lengthened  canal,  and  lengthened  canal  with  9-cu  m/sec 

supplemental  withdrawal 
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Figure  18.  Head  difference  between  the  reservoir  side  (with  9-cu  m/sec 
supplemental  withdrawal)  and  the  river  side  of  Bushy  Park  Dam 


Figure  19.  Hourly  data  for  the  intrusion  period  starting  30  May  1988 
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Figure  20.  Example  model  results  for  three  intrusion  events 


APPENDIX  A:  NUMERICAL  TIDAL  MODEL 


1.  This  appendix  presents  and  describes  the  mathematical  equations  that 
make  up  the  flow  and  salinity  portions  of  the  FIBS  model  used  in  this  study. 
Due  to  the  complicated  nature  of  this  subject,  the  description  will  be  limited 
to  an  overview.  Blumberg  (1977)*  and  Wang  (1983  and  1984)  describe  similar 
laterally  averaged  flow  and  salinity  models  based  on  the  finite  difference 
method. 


Equations 


2.  The  following  governing  laterally  averaged  dynamic  differential 
equations  describe  estuarine  flow,  mixing,  and  circulation: 

&.  The  vertically  integrated  continuity  equation. 

b.  The  horizontal  momentum  equation. 

c.  The  continuity  equation. 

The  salinity  transport  equation. 

Additional  terms  were  included  in  the  equations  to  describe  both  lateral  in¬ 
flows  and  tributaries.  The  following  algebraic  expressions  are  required  to 
dynamically  modify,  or  close,  the  set  of  differential  equations: 
fl.  The  equation  of  state, 
fe.  The  friction  coefficient  equation. 

£.  The  equations  for  vertical  eddy  diffusivities  for  mass  and 
momentum . 


An  improved  vertical  eddy  diffusivity  closure  was  made  using  a  differential 
equation  for  the  transport  of  turbulent  kinetic  energy. 

3 .  Laterally  averaged  equations  use  values  averaged  across  sections  at 
certain  increments  of  depth.  The  channel  widths  are  specified  at  each  node 
depth  as  indicated  in  the  hypothetical  channel  cross  section  (Figure  Al) . 

&■  Vertically  integrated  continuity  equation: 


3(B0h) 

8t 


*  ^  i  <ra)  “  - 


0 


(Al) 


*  References  cited  in  this  appendix  can  be  found  in  the  References  at  the  end 
of  the  main  body  of  the  report. 


Al 


Figure  Al.  Hypothetical  channel  cross  section  (looking  seaward) 
b.  Horizontal  momentum  equation: 


A (DB)  *  A (UUB)  -  BB»  IS  *  A  f™3  -  BN*  A]  <a2> 

Cd  0|U|  |»|  <“]-« 

C.  Continuity  equation  in  a  vertical  plane: 

+  *  (WB)  -  0  (A3) 


<&.  Salinity  transport  equation: 


3(SB) 

~Tt~ 


(SUB)  -  BK, 

ax2 


♦  A <SUB  - BK-  A> 


S,q 


HA,  £  B  dZ 


(A4) 
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&  Equation  of  state: 


p  *  pa  (1  ♦  C4>  •  S) 


(A5) 


where 

B„  -  breadth  at  water  surface 
h  -  deviation  in  water  surface 
t  -  time 

X  -  horizontal  distance 
R  -  depth  to  mean  tide  level 
U  -  horizontal  (laterally  averaged)  velocity 
B  -  breadth 
Z  -  vertical  distance 
q  -  lateral  inflow  rate 
A.  -  element  length 

Nx  z  -  horizontal  and  vertical  eddy  diffusivity  for  momentum 
W  -  vertical  (laterally  averaged)  velocity 
Cd  -  quadratic  friction  coefficient 
g  -  acceleration  due  to  gravity 

p0  -  freshwater  reference  density  (0.9987  g/cu  cm) 

S  -  salinity  parameter  expressed  here  as  chlorinity,  ppm 
K,  ,  -  horizontal  and  vertical  eddy  diffusivity  for  mass 
i  -  lateral  inflow  quantities 
H  -  total  depth 

-  constant  (1.3751  E-6) 


l-  Bed  friction: 


Cd(x) 


2.22  R1/3 


(A6) 


where  n  is  Manning's  coefficient  of  friction. 

Vertical  eddv  diffusivities : 

4.  Under  homogeneous  conditions,  characteristic  eddy  diffusivities  are 
NJO  -  K„0  -  K£  ,  where  Ki  is  an  specified  constant.  The  effect  of  vertical 


A3 


density  stratification  is  included  using  the  method  of  Munk  and  Anderson 
(1948).  Thus: 

K,  -  K,0  (1  +  3.33  Ri)-<3/2> 

and 

N*  -  N,0  (1  +  10  Ri)_<1/2) 

where 

Ri  *  If  .  dp/dZ 

po  (au/az)2 


Then  at  depth, 

K*(z)  -  F,(z)  •  K,  <A7) 

where  F,  is  a  similarity  distribution  for  vertical  eddy  diffusivity. 
Turbulent  kinetic  energy 

5.  Instead  of  assigning  a  value  to  KJ  ,  the  model  will  optionally  use 
turbulent  kinetic  energy  k  generated,  dissipated,  and  transported  by  the 
flow  to  calculate  vertical  diffusivities .  The  scheme  employed  is  similar  to 
that  described  by  Smith  and  Takhar  (1979)  and  Smith  (1982).  Transport  is 
computed  in  one  dimension  by 


8k 
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Pv  +  e  =  0 


(A8) 


where  Pk  is  production  of  turbulent  kinetic  energy 


Pk  =  22.4 


«  is  dissipation 


e  *1.66 
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U  is  the  cross-sectional  average  velocity 


U 


1 

5 


h 

£  BU  dz 


and 

D*  -  horizontal  diffusion  coefficient 
U*  -  shear  velocity 
A  -  cross-sectional  area 

Then,  K*0  -  Nxo  -  0.09  k1,z  R  ,  and  the  values  for  K,  ,  Nx  ,  K*(z)  and  Nx(z) 
are  computed  as  previously  described. 

Boundary  Conditions 


Ocean  boundary 

6.  In  general,  the  varying  tidal  conditions  are  controlled  at  the  ocean 
boundary  of  the  model  mesh.  This  is  accomplished  by  specifying  a  tidal  water 
level  target  value  at  each  time-step  in  the  form  of  h'  (t)  |  x_0  -  /( t)  .  The 
model  has  the  capability  to  synthesize  tidal  sequences  for  any  specified 
starting  time  using  12  constituents  and  harmonic  coefficients  provided  by  the 
National  Ocean  Survey.  The  individual  constituents  used  included  M2,  S2,  Kl, 
01,  PI,  N2 ,  L2 ,  Sa,  M,f,  Mf,  Mm,  and  Ssa. 

7 .  A  nonreflecting  ocean  boundary  was  computed  for  this  study  by  re¬ 
placing  the  vertically  integrated  continuity  equation  with  the  following 
expression: 


_  |h(t)  -  h'  (t)  j 


(A9) 


where 


Cw  -  (gH)1/2|x.c 
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and  T£  is  a  damping  parameter  (found  to  be  about  500  sec).  This  method, 
similar  to  that  presented  by  Blumberg  and  Kantha  (1985),  allows  long  waves,  or 
water  level  perturbations,  originating  within  the  model  domain  to  pass  through 
the  ocean  boundary  without  reflection.  The  nonreflecting  ocean  boundary  per¬ 
mits  faster  stability  within  the  model.  In  this  study,  the  model  reached  the 
stable  condition  in  two  tidal  cycles.  During  the  flood  tidal  phases,  salinity 
concentration  was  specified  at  the  ocean  boundary;  otherwise  the  transport 
boundary  was  unconstrained. 

Upstream  boundary 

8.  Velocities  are  specified  at  upstream  inflow  boundaries  and  flow 
rates  are  specified  at  lateral  inflows.  In  this  study,  the  condition 
U(t)  -  0  was  specified  at  Pinopolis,  the  upstream  end  of  the  model,  and  a 
lateral  inflow  was  specified  at  the  adjacent  element  to  represent  the  hydro- 
power  flow  release.  The  Cooper  River  contributes  almost  all  of  the  freshwater 
inflow  to  the  system  and  is  controlled  at  the  Pinopolis  Dam.  Salinity  concen¬ 
trations  (usually  zero)  were  specified  at  these  inflows. 

Bottom  boundary 

9.  For  this  study,  the  bottom  boundary  (Z  -  0.0)  was  considered  to  be  a 
slip-flow  boundary.  The  shear  stress  was  imposed  by  the  equation; 


N,  ™  =  Cd  UJU| 


(A10) 


Optionally,  the  model  will  extrapolate  shear  stress  linearly  from  within  the 
flow  to  the  bottom  boundary  using  a  method  suggested  by  Smith  (1982).  Within 
the  salinity  transport  equation,  no  flux  is  allowed  at  the  bottom  boundary  and 


us  ♦  *■  S  - 0 


(All) 


Vertically  Stretching  Coordinates 

10.  In  the  vertical  or  Z-direction,  computational  nodes  are  assigned  to 
the  surface ,  intermediate ,  and  bottom  levels .  When  changes  occur  in  the  water 
levels,  all  nodes  except  the  bottom  node  move,  or  stretch,  to  conform  with  the 
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depth  of  flow.  The  typical  configuration  of  the  nodes  within  an  element  are 
shown  in  Figure  A2.  The  X  and  Z  are  the  actual  world  coordinates.  The  x 
and  z  are  the  nondimens ional  element  coordinates  that  vary  in  each  element 
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Figure  A2.  Example  computational  ele¬ 
ment  from  schematic  mesh  showing 
stretching  coordinates  and  nodes 


from  0.0  to  1.0  upstream  and  vertically  upward,  respectively.  The  two 
coordinate  systems  are  related  by 

8  18 

\ 32  H  Hz 


a  m  l  a 
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(A12) 

(A13) 


The  stretching  coordinates  system  introduces  new  terms  into  the  governing 
differential  equations.  A  derivative  in  some  level  of  the  X-direction  is 
calculated  by  the  following  expression: 


a  a  _  [  z  ah  a  +  (i  -  z)  aRl  a 
EX  A.ax  Ira,  Ex  EZ  Ra,  ExJ  EZ 


(A14) 
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This  method  is  similar  to  the  finite-difference  methods  presented  by  Sheng 
(1983). 


Solutions  of  the  Equations 

11.  The  governing  equations  were  solved  over  the  previously  described 
finite  elements  using  an  explicit  orthogonal  collocation  method  of  weighted 
residuals.  Weighted  residual  methods  can  be  used  to  expand  an  unknown  solu¬ 
tion  in  a  set  of  trial  functions  T  over  a  domain  v  .  The  collocation 
method  applies  weighting  functions  (wj  ,  a  Dirac  delta  function)  to  the  trial 
functions ,  chosen  to  be  a  set  of  orthogonal  polynomials .  The  node  locations 
Xj  satisfy  the  roots  of  the  polynomials.  The  form  of  the  collocation  method 
is  now 


£  WjT  dv  * 


(A15) 


where  Wj  -  6(x  -  Xj)  . 

12.  Instead  of  solving  for  the  unknown  trial  function  coefficients,  a 
major  simplification  described  by  Finlayson  (1972)  is  that  the  solution  can  be 
derived  in  terms  of  values  at  the  collocation  nodes  (r) .  The  nodes  are 
located  at  Gauss  points  for  this  purpose ,  and  the  method  is  known  as  ortho¬ 
gonal  collocation.  The  derivatives  and  integrals  can  therefore  be  expressed 
in  terms  of  values  of  the  function  at  the  collocation  nodes  by 


J»Ar  and  ■  B  r  (A16) 

OX  dx2 

The  matrices  A  and  B  are  computed  from  the  trial  functions. 

Integration  in  Time 

13.  An  explicit  time  integration  method  was  used  to  compute  variables 
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dynamically  over  time.  The  fourth-order  Runge-Kutta  method  was  used  in  this 
study.  This  method  is  very  accurate,  and  often  used  to  solve  "stiff"  equation 
sets,  equations  with  distinct  stability  limits.  Results  were  found  to  be 
independent  of  time-step.  At,  up  to  a  Courant  number,  At/A.-Cw  ,  of  1.0  (or 
2.0  based  on  longitudinal  node  spacing).  Starting  with  the  initial  conditions 
yn  at  time  t„  ,  the  function  y'  (t)  -  f(t,y)  where  y'  is  the  rate  of 
change  in  y  is  solved  in  four  steps  over  At  : 

a.  Kq  *  At  •  f  (tn.yj  . 

b.  K4  -  At  •  f  ^  ^  ,  yn  *  . 

c.  K2  -  At  ■  f  jt„  ♦  ^  ,  yn  ♦  . 

d.  Kj  ■  At  ■  f  |t„  *  At  ,  yn  +  K2j  . 

Finally,  the  result  of  the  four  steps  are  combined  to  complete  the  time-step: 

y»*i  -  yn  ♦  l  (Ko  ♦  2k,  ♦  2K2  ♦  K3)  (A17) 
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